Selection theory of free dendritic growth in a potential flow 
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The Kruskal-Segur approach to selection theory in difFusion-limited or Laplacian growth is ex- 
tended via combination with the Zauderer decomposition sclieme. This way nonlinear bulk equations 
become tractable. To demonstrate the method, we apply it to two-dimensional crystal growth in 
a potential flow. We omit the simplifying approximations used in a preliminary calculation for the 
same system [T. Fischaleck, K. Kassner, EPL 81, 54004 (2008)], thus exhibiting the capability of 
the method to extend mathematical rigor to more complex problems than hitherto accessible. 
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I. INTRODUCTION 

Pattern formation is ubiquitous in nature. Snowflakes 
constitute an everyday-paradigm of a self-organized 
structure, apparently the first that was the subject of 
scientific study [ij. Any physical pattern possesses at 
least one characteristic length scale, and if it is dynamic, 
it also has a characteristic time scale. The foremost task 
of scientific endeavour in the field of pattern formation is 
to explain the emergence of these scales and to determine 
them quantitatively. Since systems with linear dynam- 
ics will, due to the superposition principle, not normally 
single out a particular length scale, an essential ingredi- 
ent of pattern- forming systems is the nonlincarity of their 
dynamics [1]. 

As it turns out, snowflakc-like structures - dendritic 
morphologies - also arise at microscopic scales in the 
casting of metals and they determine structural proper- 
ties such as the strength of the material, which imparted 
considerably more importance to scientific preoccupation 
with them than just fundamental interest would have. 

The first models of dendritic crystal growth assumed 
transport of heat away from or material to, the grow- 
ing nucleus to be simply diffusive, so the term difFusion- 
limited growth was coined. Understanding the selection 
of dynamical features such as a basic length scale and 
the growth velocity turned out to be remarkably difficult 
even within these simplifying models. Almost forty years 
passed between Ivantsov's approximate solution (3| that 
did not exhibit selection and the development of an ana- 
lytic theory explaining the mechanism of structure selec- 
tion • This may seem even more surprising consider- 
ing that the bulk equations of difFusion-limited systems 
are linear and the nonlincarity of the dynamics emerges 
solely via the equations of motion for the two-phase in- 
terface. In fact, the analytic approaches developed had 
to rely heavily on this linearity. 



[Klaus. Kassner@ovgu.de| 



Ivantsov's theory, neglecting surface tension at the 
boundary between the melt and the solid, predicts only 
the product of the tip radius of a dendrite and its growth 
velocity, as a function of the undercooling. In experi- 
ments, the undercooling determines both quantities sep- 
arately. A decisive step towards the solution was the 
insight that without surface tension the problem is ill- 
posed 0, Q and that the capillary length has to be taken 
into account, even if it is much smaller than any length 
scale of the arising pattern. Surface tension regularizes 
the mathematical problem and drastically alters the solu- 
tion space. Without surface tension, there is a continuum 
of parabolic needle crystal solutions. With isotropic sur- 
face tension, there are no solutions with a shape close 
to one of these Ivantsov parabolas (or paraboloids), no 
matter how small the surface tension, a fact that testifies 
to the singular nature of the "perturbation" surface ten- 
sion. With anisotropic surface tension, the continuum 
oF Ivantsov solutions is reduced to a discrete set, with 
the Fastest oF the needle crystal solutions being the only 
linearly stable one @. Hence, the selection problem is 
broken down into two parts - an existence problem for 
a discrete set of solutions and the stability analysis sin- 
gling out one element of the set as the one that should be 
observed. A completely analogous theory was developed 
For Saffman- Taylor fingers in viscous fingering |10l4l2| . 
where selection is also due to surface tension, albeit not, 
of course, to its anisotropy. 

These theories were two-dimensional, just as the origi- 
nal numerical work giving evidence for a selection mech- 
anism based on solvability [l3l - [l5| . Three-dimensional 
situations considered initially referred to axisymmetric 
crystals [H [3, hence were not very realistic. Later, 
steps were taken to extend the theory towards non- 
axisymmetric needle- crystal shapes ^vA [isj and eventu- 
ally, an analytic th eory was developed for the fully non- 
axisymmetric case [l9l [20j ; all the necessary elements of 
the final complete theory were not present before ReF. [2O. 

From the outset, two different analytic approaches 
were pursued. With the first method, the equation of mo- 
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tioii of the two-phase interface is Hnearized, which leads 
to an intcgro-difFerential equation in non-local problems 
(such as dendritic growth or viscous fingering). Using 
Fredholm's alternative, a solvability condition is derived 
that is satisfied only by a discrete set of values of the se- 
lection parameter (a nondimensionalized surface tension 
or its inverse) [2ll |. The second approach, pio- 

neered by Kruskal and Segur [22| , consists in solving the 
interface equation far from singular points in the complex 
plane via a perturbation expansion in terms of the small 
selection parameter and in the vicinity of these points 
via a scale transformation and reduction to a local equa- 
tion. The two solutions then have to be asymptotically 
matched to obtain a globally valid solution. Parame- 
ter relationships established in the matching procedure 
yield the selection criterion [23l - l26j . There is general 
agreement that only the second approach is mathemati- 
cally rigorous [l^, HI, 113, HI] ■ The linearization of the 
first method introduces approximations that normally 
will not invalidate the scaling relations obtained; but it 
will not reproduce their prefactors correctly nor provide 
a globally valid approximate solution. Moreover, if the 
equations contain more than one small parameter (say, 
a kinetic coefficient or a characteristic number describ- 
ing the flow, besides the usual selection parameter), the 
linearization may produce even worse results due to the 
structural instability of the problem (23j . 

In both approaches, it is necessary to first derive an 
(integro-differential) equation for the interface position 
depending on a single independent variable. This can be 
achieved, e.g., by eliminating the bulk field variables via 
conformal mapping (in the viscous fingering case) or us- 
ing Green's function methods (in crystal growth). These 
techniques are applicable only for linear bulk equations, 
which seemed to preclude utilization of the method for 
convection problems. 

For a long time, the only exception to this restriction 
has been the work on two-dimensional crystal growth in 
an Oseen flow by Bouissou and Pelce [29[. To obtain 
the selection criterion, they used a method, outlined in 
Ref. HO, which is closely related to the first of the two ap- 
proaches mentioned, hence not rigorous. Their method 
has recently been extended by Alexandrov et al. [IH to 
include solute diffusion. The equation of motion for the 
deviation of the solution from that of the problem with- 
out surface tension is simplified in the style of a linear sta- 
bility analysis which allows to avoid the derivation of an 
integro-differential equation. Moreover, the adjoint lin- 
ear operator is constructed hcuristically in Fourier space 
to obtain a solvability condition in the spirit of the Fred- 
holm alternative, a procedure that may introduce addi- 
tional (possibly problematic) approximations. 

A different method, having the potential of achieving 
the same level of rigor for problems with nonlinear bulk 
equations as the asymptotic matching approach, was re- 
cently introduced j32l |33| . It consists in a combination of 
Zauderer's decomposition scheme [s^l for partial differ- 
ential equations with the Kruskal-Segur approach. Zaud- 
erer decomposition is the step allowing reduction of non- 



linear bulk equations to an interface equation and thus 
circumventing the necessity of an exact integral equation, 
available only for linear bulk equations. Reference |32| . 
dealing with potential flow, was more or less a proof of 
concept, in which we copiously used additional approx- 
imations to simplify the result to an easily digcstable 
form, allowing to map it in the end to the fiowless fi- 
nite Peclet number case treated by Ben Amar [23|. The 
main purpose of the present paper is to remove these 
approximations, which renders the treatment more com- 
plex, but does not impose insurmountable obstacles. Of 
course, the mapping obtained gets lost, because it was 
only approximate. Clearly, potential flow is not a very 
realistic assumption, but it has the advantage that the 
unperturbed problem is exactly solvable; the analog of 
Ivantsov's analytic solution exists. This is different in the 
case treated by Bouissou and Pelce [1^, where already 
the zeroth-order problem is solved within an approxima- 
tion, replacing the Navier-Stokes equations with the Os- 
een problem (which in two dimensions does not even give 
a uniform approximation to the true flow [HI). We will 
consider more realistic flows and a better approach than 
the Oseen approximation in a later publication. 

The paper is organized as follows. In SecHIl the model 
equations are given. They are nondimensionalized and 
rewritten in parabolic coordinates in Sec. IIIII Section HVl 
gives the analog of the Ivantsov solution in the presence 
of a potential flow. Then the method of Zauderer de- 
composition is explained in Sec. |Vl allowing us to reduce 
the set of partial differential equations of the full problem 
to an integro-differential equation for the interface alone, 
without losing the terms decisive for solvability theory. 
Next, the decomposed equations are solved to flrst order 
(in the small parameter cr^^^) in Sec. I VII Near the solid- 
liquid interface, the relevant behavior beyond all orders 
of regular perturbation theory is obtained in Sec. IVIII 
from a WKB analysis. On the other hand, in Sec. IVIIIl 
the asymptotic Kruskal-Scgur reduction to a locally valid 
equation, applicable near a singularity in the complex 
plane, is carried out. It leads to a nonlinear integro- 
differential equation constituting an eigenvalue problem. 
The numerical solution of this eigenvalue problem deter- 
mines the selected velocity (and other properties) of the 
needle crystal. Detailed results are given for a set of pa- 
rameters corresponding to a particular experimental sys- 
tem, which however does not exhibit potential flow, so 
the comparison is only qualitative. Some conclusions are 
offered in Sec. IIXI A few general calculations and slightly 
elaborate mathematical conversions are relegated to two 
appendices. 

II. MODEL EQUATIONS 

Heat transport in the liquid and solid phases is 
described by the diffusion-advection equations 

rIT 

^ + {^.V)T = DV'T (1) 
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with w the flow velocity in the liquid and w e in the 
solid. The advection term coupling the temperature and 
flow equations renders the bulk problem nonlinear, de- 
spite the simplifying assumption of potential flow made 
below (which reduces the flow description to a linear 
equation). As the notation suggests, we assume the ther- 
mal diffusivity D to be the same in both phases (symmet- 
ric model). We consider an incompressible flow, which 
means that a stream function can be introduced. In two 
dimensions, its defining equation takes the form 



(2) 



Because depends on x and y only, we obtain V x w 
-V^V'G^ and taking the flow to be potential, we have 



0. (3) 

When specializing Eq. ([1]) to one of the phases, we 
will denote the temperature variable by T' and T'* in the 
liquid and the solid, respectively. Equations ([2]) and ([3]) 
are the bulk equations of motion for the flow. 

Because we are looking for steady-state solutions, we 
need not prescribe detailed initial conditions. We must 
however specify boundary conditions for each of the bulk 
equations. At infinity in either the liquid or solid we re- 
quire homogeneous Dirichlet boundary conditions for the 
temperature fields, i.e., we set the temperature constant: 



T'(x) 

r^(x) 



Tm 



for d(x,r) 
for d(x,r) 



(4a) 
(4b) 



Herein, F is the interface and (i(x, F) denotes a signed 
distance function, increasing towards the liquid. Tm is 
the bulk melting temperature of the solid, and for crystal 
growth to occur in a pure system, we must have Too < Tm- 
The dimensionless parameter characterizing this under- 
cooling is 



Tm - To, 
L/cp 



(5) 



where L and Cp are the latent heat and specific heat, both 
referred to a unit volume. 

Moreover, there are boundary conditions at the inter- 
face, reading 



L 



T \r = Tm datia{9) , 

Cp 

LVn = Dcp[vT\-VT\]-n. 



(6a) 
(6b) 

(6c) 

Equation (|6aP describes continuity of the temperature 
across the interface, Eq. (j6bp is the Gibbs-Thomson con- 
dition giving the equilibrium temperature of a melt- 
crystal interface with curvature k. That is, we assume 
kinetic effects to be negligible, implying local thermal 
equilibrium at the interface. If the interface is given by 

y = h{x), then k = -h" {x) I {l + h' {x)'^)^^^ . do is the 
average capillary length, 

, TmCp , , 

do =70-^, (7) 



where 70 is the angular average of the orientation depen- 
dent surface tension 7(6'). 6 is the angle of the inter- 
face normal with some fixed direction, for example the 
direction of the y axis. From the thermodynamics of 
interfaces we know [36j that the angular dependence of 
the non-averaged capillary length is not that of the sur- 
face tension but that of the surface stiffness 7(6*) +j"{9). 
This is described by the factor a{0) = {-f{9) + j"{9))/jo- 
For simplicity, we will assume fourfold anisotropy here, 
described by a single harmonic, i.e. 



a{e) = l-l3coa{A9) . 



(8) 



From a mathematical point of view, Eqs. (|6al6bp to- 
gether with the boundary conditions at infinity and some 
initial condition for the two temperature fields are suf- 
ficient to solve the diffusion problem in the two phases, 
with given flow field and interface position. However, 
the interface position is a priori unknown, its determi- 
nation is part of the problem. Therefore, an additional 
interface equation is needed. This is the third boundary 
condition, the Stefan condition ([Sc)) . Physically, it fol- 
lows from energy conservation across the interface. Vn 
is the interface normal velocity and the left-hand term 
describes latent heat production due to advancement of 
the interface, whereas the right-hand side gives the sum 
of the heat currents into the solid and liquid phases. The 
normal vector n points from the solid into the liquid. 

The flow field w is dynamic only in the liquid. We 
need a boundary condition at infinity, where we impose 
a constant flow directed opposite to the growth direction 
of the needle crystal: 



w(x) 



for d(x, F) ^ 00 



(9) 



Because we have potential flow, we cannot impose condi- 
tions for all three components of the flow velocity at the 
interface; the Laplace equation for the stream function 
does not admit prescription of more than one scalar quan- 
tity. Physically, this is reasonable, since potential flow is 
frictionless, hence we cannot prescribe the tangential ve- 
locity at the interface, there is no no-slip condition. The 
normal velocity, on the other hand, follows from mass 
conservation. We assume the simplest case, viz. equal 
mass densities in the solid and the liquid. Then the liq- 
uid is neither sucked towards the solid (which would be 
the case if the density of the solid were higher than that 
of the liquid), nor ejected from it. Since the solid does 
not move [33| , the normal velocity of the liquid must be 
zero: 



n • w(x)|r = . 



(10) 



Equations ([T]) through (fTO]) . supplemented by initial 
conditions for all the fields, constitute the complete 
mathematical description of an idealized physical system. 
Requiring the solution to be stationary and to correspond 
to a crystal growing at constant velocity V along the y 



direction we may replace Vn in Eq. ([6c|) by Vey ■ n = Vn^ 
Transforming to a moving frame of reference. 



r ^ r + Vte, 



V ' 



(11a) 
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w + Ve 



■y 1 



(lib) 



in which the interface is at rest, all time derivatives dt get 
replaced by -Vdy. Note that due to the transformation 
of w, Eq. ^ is invariant under this change of frame. 
Nevertheless, the time derivative can be dropped after 
the transformation, because we seek a time-independent 
solution. Moreover, the "flow velocity in the solid" be- 
comes -Vey by virtue of the transformation. 



III. PARABOLIC COORDINATES AND 
NONDIMENSIONALIZATION 

A family of exact analytic steady-state solutions to 
the model equations exists for vanishing capillary length 
and, similar to Ivantsov's solution in the flowless case, 
the crystal interface is parabolic [ssl - lioj . Therefore, it 
is useful to introduce parabolic coordinates. We employ 
conformal parabolic coordinates 

1 



x = v^, y = ^{v-e), (12) 

their advantage being equality of the scale factors = 

9n - ''f' ■ appendix, some of the transfor- 

mation formulas and a graphical representation of the 
coordinate lines are given. 

To nondimensionalize the equations, we use the tip ra- 
dius p of the Ivantsov like solution, defined as the inverse 
of the curvature at the tip, as a length scale. The corre- 
sponding diffusion time p^/D is taken as a time scale 

K p^ 
x,y ^ px,py, K-> - , t^—t. (13) 
P D 

Note that this implies ^ and 77 to scale with ^/p: £,,r] ^ 

The nondimensional form of the flow velocity follows 
immediately 



D 

w ->• — w . 

P 



(14) 



(This implies ip Dip.) The nondimensional flow veloc- 
ity at infinity then becomes 



Pf 



pU 
D 



(15) 



the so-called flow Peclet number. Moreover, it will turn 
out useful to include the growth Peclet number 



Pr 



pV 
D 



(16) 



into the prescription for nondimensionalization of tem- 
perature 



L 



(17) 



which means that the nondimensional temperature in the 
liquid approaches -A/Pc at infinity. 

With these transformations and dropping time deriva- 
tives, as we are interested in stationary solutions, we ob- 
tain the following set of bulk equations 

^^tJ-^^T^=tJ^ + t4,, (18a) 
Pc{^T^-vT^)=T^i + T^r,, (18b) 
V'?? + i^vv = . (18c) 

The boundary conditions for the fields at the interface 
now read 



T' = --<JKa{e) , 



(19a) 
(19b) 



[CVs]' = -v's{Ti-T^)+T^-T^, (19c) 

i^i + V'.^v = Pc iVs + v'A) , (19d) 

where ?7s(^) is the interface position, a prime means a 
derivative with respect to ^ along the interface, and 



2do 

pPc 



(20) 



is the selection or stability parameter. (In order not to 
overburden the notation, we have dropped the qualifier 
|r next to the fields and their derivatives, indicating that 
these quantities have to be evaluated at the interface po- 
sition.) 

Finally, the boundary conditions at infinity may be 
written (see Fig. O 



(7? 00) 

Pc W ; 



(21a) 

- (7? < 1, 1^1 ^ 00) , (21b) 

i>-{Pf + Pc)ri( iv^^), (21c) 

the last equation being an asymptotic equality (4T| . 



IV. EXACT SOLUTION IN THE ABSENCE OF 
SURFACE TENSION 

If wc neglect surface tension in the spirit of Ivantsov, 
eqs. p9all9bp tell us that with vanishing capillary length 
the interface becomes an isotherm, T'|r = T*|r = 0. As- 
suming it to be a coordinate line suggests the tempera- 
ture field to depend on one of the coordinates only, i.e. 
T = T{ri). Inserting this into (|18all8bp . we have 



-Pcvt:, 



rjil 



(22a) 
(22b) 



In view of boundary condition (|21bp . wc see that the 
second equation is solved by = 0. The first can have 
a purely rj dependent solution only if 



(23) 
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Inserting this into (|18cp , we find 

f"{n)=Q /(??)= Cl77 + C2 



(24) 



with constants ci and C2 to be determined from the two 
boundary conditions (jl9dl21c|) . By assumption, the in- 
terface is at (^) E 1, so (|19dp simphfies into V'c(l) = Pc- 
We find ci = Pf + and C2 = -P/. Therefore, 



^ = v'^ = e(Pc?y + P/(?7-i)). 



(25) 



This leaves us with the ordinary second-order differen- 
tial equation ((22a|) T^^ + {P^-q + Pf{ri - 1)) T/, = subject 
to the boundary conditions T^l) = and (|21ap . The 
solution is straightforward: 



drj', (26) 



where Pc is determined as a function of A and Pf from 



P^ 



/ 



-^(v-i) 



dr/'. 



(27) 



In the limit Pf ^ 0, this becomes identical to the usual 
Ivantsov relation for diffusion-limited dendritic growth, 
whereas for Pc « Pf we can evaluate the formula an- 
alytically, which yields A = Pc^TT/2Pf . The selection 
problem arises in the same way here as in the flowless 
case: only Pc is determined, given the undercooling and 
the imposed flow, by the Ivantsov-like solution. Given 
Pc, we can calculate the product of the growth rate V 
and the tip radius p but not both quantities separately. 
Ivantsov's approach eliminates the parameter a from the 
equations. As (|20p tells us, we may calculate p once we 
know a and Pc- So the aim must be to include a into the 
theory and to obtain a value for it. 

For easy reference, we will call the solution (|25p - fI7} 
the flow-Ivantsov solution. 



V. ZAUDERER DECOMPOSITION AND 
CONTINUATION TO THE COMPLEX PLANE 

The approach to be followed is singula.r perturbation 
theory about the flow-Ivantsov solution [43 • Thus, we 
shall consider small deviations from it 



and V V'"'"^ + "0 J 



(28) 



but we will be careful to avoid illegitimate lineariza- 
tions. In particular, we will not linearize terms containing 
derivatives of the interface position. 

In this first approach, we shall restrict ourselves to the 
limit of small growth Peclet number, i.e., Pc « 1- As it 
turns out, interesting results for the selection parameter 
arise only, if we also assume Pc « Pf. In particular, this 
means that terms with a factor of Pc will be neglected in 
the exponentials of Eqs. ([^S]) and (P7|) . 



The diffusion-advection equation (|18ap then becomes 
an inhomogeneous equation for the temperature devia- 
tion, in the liquid, from the flow-Ivantsov solution 

rJ^ + T^,-(V„ + ^P/)Ti 

+ (^e + P/(7?-l))7;; = ^^e-*(''-i)\ (29) 

The only approximation in this equation is that Pc has 
been set equal to zero. All fleld nonlinearities are still 
present. 

In the solid (Eq. (jl8b[) ). we just drop the terms multi- 
plied by Pc and obtain a Laplace equation 



7|« + r;, = 0, 



(30) 



and Eq. pScp . being linear, remains formally unchanged. 



+ V'r, 



0. 



(31) 



but the meaning of ip is now different (it is the deviation 
of the stream function from its form in the flow-Ivantsov 
solution) . 

To rewrite the boundary conditions, we set Ty^ = 1 + 
/i(^). While it is legitimate to view /i(^) as a small quan- 
tity, this may not be true for h'{^) and higher deriva- 
tives. Expanding (jl9ap about the flow-Ivantsov solution, 
we have T^^\r + T^\r = d^T^\^ih + T'(^,l) = T\ = 
T*(^, 1). Later, we will need the derivative of the in- 
terface temperature with respect to ^. We are then 
not allowed to simply take the partial derivative of the 
temperature field with respect to which would give 

- h'{C) = Tl (because dr,T^^\r,=i = -1). Actually, the 
derivative must be taken along the interface, so we obtain 
T\ + T,\h'{0 - h'{^) = + r,^/i'(0 instead. Keeping this 
proviso in mind, we may use the simpler form before dif- 
ferentiation as the shortest description of the appropriate 
boundary condition. The full set of interface boundary 
conditions then reads 



T" + h = T , 



Ip^ + h'ipr, 



-Pfm' 



(32a) 
(32b) 

(32c) 
(32d) 



All of these are evaluated at 77 = 1, but if they are to 
be differentiated, then this has to be done before setting 
77 = 1 and field derivatives with respect t o ri will arise. 

Zauderer's asymptotic decomposition [3J] is a projec- 
tion scheme reducing the solution of a system of partial 
differential equations (PDEs) to the solution of a series 
of first-order equations. These may be decoupled within 
a perturbative approach, if a small parameter or slow 
variable is available. (Otherwise, the series of first-order 
equations does not offer any simplification over the origi- 
nal system of PDEs.) Originally conceived for hyperbolic 
equations, the method generalizes to elliptic systems in 
the complex plane. 
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Asymptotic decomposition seems to have largely 
passed into oblivion (at least in the physics world) , possi- 
bly because often multi-scale expansions are superior to 
it, leading to more easily tractable equations. Neverthe- 
less, for the problem considered here, Zauderer decom- 
position is particularly well-suited, not losing informa- 
tion about transccndentally small terms (if the "principal 
part" [34 of the set of equations is chosen correctly). 

In [3^ , we discussed that in the case of purely diffusive 
transport and in the limit of small growth Peclet num- 
ber, Zauderer decomposition of the transport equations 
is equivalent to their factorization 

{d^ + id^)T' =0, {d^-id^)T' =0, (33) 

and that a few simple manipulations of these partial dif- 
ferential equations using the boundary conditions lead to 
a local equation for the interface position /i(^), 

aa{e)K = {l-iOh{0. (34) 

This equation contains the complete information needed 
to compute the mismatch function that has to be zero at 
the tip of the needle crystal for selection to be possible. 
Near the singularity ^ = -i in the complex plane, Eq. ([34]) 
describes the dominant behavior of the solution. 

Equations (|33|) are, due to their simplicity, well-suited 
for a discussion of the strategy of our approach. Typ- 
ically, Zauderer decompostion produces, from the basic 
partial differential equations of the problem, a leading- 
order or "principal-part" equation in each domain that is 
first order and easily solvable (for example by the method 
of characteristics). The complete solution will be a sum 
of this leading term and other contributions that may 
be calculated in subsequent steps of the perturbative 
scheme. In the simple case considered here for explana- 
tory purposes, the temperature field satisfies the Laplace 
equation factorizing in the complex plane (a formal Za- 
uderer decomposition just reproduces this factorization). 
Solving the factorized equations, one finds 

m, V) = hii + Kv - 1)) + /2(C - - 1)) (35) 

with analytic functions /i and /2. Inserting this solution 
into the boundary conditions at 77 = 1 and analytically 
continuing the resulting equations to w = -i (where w is 
the analytic continuation of ^), some of their terms must 
become singular in the limit h{Q « 1 (ct ^ 0), because 
the curvature term in ([34]) becomes singular in that limit 
[see Eq. (jA9[) ]. The solution for T' must be analytic in 
the liquid, i.e., for > 1, corresponding to the upper half 
complex plane. Hence the /2 term will remain analytic 
near w = -i [43j . The important contribution to T' that 
may diverge near the singularity will then come from /i , 
and this function is the solution to the first equation of 
(|33)) . Similarly, it may be argued that in the solid the /2 
term is the important one; it solves the second equation 
of ([55]) . Thus, Eqs. ([55]) give us a solution that is valid 
near the singularity. At the interface, far away from the 
singularity, this approximation is also justified (albeit to 
a lesser degree), because there the curvature term in the 



boundary conditions can be linearized. For the (homoge- 
nous part of) the linearized equations, each of the /i and 
/2 terms alone is a solution. To obtain the full solution, 
the second singularity at w = i has to be taken into ac- 
count. To treat that case, the lowest order equations for 
r' and T'* would have to be interchanged. So the choice 
of the sign of the 9,, term in the factorization depends 
on the singularity considered, and the signs for the liquid 
and solid domains must be opposite to each other. In 
our simple example, it is sufficient to just consider one 
singularity, because the result for the second will be the 
complex conjugate of that for the first. 

Having an asymptotic solution that is valid both near 
the singularity and all the way to the interface, we may 
then impose the solvability condition of vanishing mis- 
match function on this solution j4^ . 

In the general case, we cannot simply factorize the ba- 
sic partial differential equation but achieve the reduction 
of order enabling analytic solutions by Zauderer decom- 
position. Analytic continuation to the complex plane will 
again prove useful. A convenient starting point consists 
in representing the partial differential equations cLS 0j set 
of first order equations. We define 




UO = + ^Pf , Vo = -tfj^ - Pf {ij - 1) , 

Fo = -e-^^"-!)' . (36) 

The governing equations then become 

W{ + AW^ + BW + CV = (37a) 
W| + AW^ = (37b) 
+ AVr, = (37c) 

A is a constant matrix, B and C are assumed to vary 
slowly as functions of rj and ^ in the vicinity of the 
Kruskal-Segur point {w = -i). This suggests a scale trans- 
formation ^, r] ^ eS., sr], emp hasizing the derivative terms 
in (|37p . As discussed in [321, ^ ^"6 related to the sta- 
bility parameter after solution of the selection problem, 
giving s oc a^l'^ . We will expand equations in powers of 
e, drop terms of order and higher and set e = 1 af- 
terwards to simplify the notation. A key of Zauderer's 
approach is to rewrite the system of equations in terms 
of eigenvectors of the matrix A appearing in its princi- 
pal part (here given by expressions of the form + Af^) 
and to obtain decoupled equations for their coefficients, 
using appropriate projections onto the eigenvectors. The 
eigenvectors of A are 

rr^(-;), r2^(;), (38) 

corresponding to the eigenvalues i and -i, respectively. 
Since A is antihermitean and the eigenvalues different. 
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these eigenvectors are orthogonal (but the formahsni does 
not rely on this). We expand the field vectors in terms 
of ri and r2 



W = Mri + eNY2 
W = N''V2 , 
V = xri , 



(39a) 
(39b) 
(39c) 



where the choice of a prcfactor s in front of the coefficient 
function N is dictated by our expectation of this term be- 
ing small in the liquid, because ri is the eigenvector lead- 
ing to an equation of the form + iM^ = in the limit 
e that can be identified with the flowless case (see 
Eq. (|37ap . where the B and C terms become negligible 
after the scale transformation in the limit e ^ 0). For the 
case without flow, we have identified this form to corre- 
spond to the equation generating the relevant component 
of our solution in the vicinity of the singularity w = -i. 
With fiow, there will be corrections of order e that we 
wish to calculate. That we have completely dropped one 
of the eigenvectors in Eqs. (|39bp and p9cp is due to the 
fact that Eqs. (j37bp and p7cp have only principal parts, 
so the coefficients of the dropped eigenvectors decouple 
completely. Equation (|37b[) holds in the solid, so we ex- 
pect the r2 term to be dominant, Eq. (|37cp refers to the 
liquid domain, so the ri term should be dominant. Once 
we have calculated the coefficient functions M, N, A^", 
and Xi may obtain the temperature and flow fields 
from 



-i{M-eN) 



rpS 



M + eN , 
X, 



(40a) 
(40b) 
(40c) 



equations that also allow us to obtain boundary condi- 
tions for the coefficient functions from Eqs. ([5^ . 

Plugging Eq. (|39|l into Eq. ((37|) and neglecting terms 
of order e^, we find 



M^ri + eN^r2 + iM,,ri - ieA^'j^ra 
+eBMri + eCx^i 



(41a) 
(41b) 
(41c) 



X? + iXr, 

Next, we project these equations onto the eigenvectors 
to cast them in the simplest possible scalar form. Pro- 
jection operators on the two eigenvectors are easily con- 
structed by tensorial multiplication with the dual vectors 
of the biorthogonal system constructed from ri and Y2- 
We have 



P2 

P^rk 



2 



(42) 



= 1,2 . 



Abbreviating a = \ (mo + ii'o), we may write the projec- 
tions of BYi and Cr,; as 



PiBvi = -ari , PiCvi = — ri , 
P2Bri = ar2, P2Cri = -— r2 



(43) 



Applying the projection operators to (|41ap and setting 
£ = 1, we obtain 



+ iM,, - aM +^x = 0, 
N^-iN^ + aM-^X = 0, 



(44a) 
(44b) 



which together with (j41bp and (|41cp gives us four equa- 
tions for the four quantities M, N, Ns, and x- 

Finally, the boundary conditions ([32]) have to be trans- 
formed into boundary conditions for our new fields. Be- 
cause W, W", and V are defined in terms of derivatives 
of the temperature field and the stream function, we have 
to take derivatives in (j41c[) . wherever non-differentiated 
temperatures appear. It is here, where care has to be 
taken that the boundary conditions hold along the in- 
terface and hence we do not obtain a boundary con- 
dition for directly from p2ap or (j32bp but one for 
dT/dC = Tj + /I'T^. Using (gOl), we find the interface con- 
ditions: 



M 

N-N' 



2 1 + ih' ' 
2 1-ih' 



_ i_a{Ka{0)y 
2 1-ih' 



X ■■ 



i+ih' 



(45a) 
(45b) 
(45c) 
(45d) 



where the prime always denotes a derivative with respect 
to ^. Combining the second and third equations, we see 
that the equations for M and N decouple from that for 
A^", because we can give their boundary condition at the 
interface without solving the equation for explicitly. 
(Of course, we have to make sure that there is a solution 
in the solid, so the behavior of N'^ near the singularity 
in question is important.) 

The boundary conditions at infinity follow from (j2ip 
combined with (j28p and simply require all fields to go to 
zero in the appropriate infinite limit: 



A/^0 
N ^0 

X-0 



{rj oo) , (46a) 

(?7 ^ oo) , (46b) 

iv < 1, lei - °°) , (46c) 

(?7 ^ oo) . (46d) 



The alert reader may be surprised that we have eight 
boundary conditions [Eqs. (|45p and (|46|) ] for four first- 
order difl:erential equations [Eqs. (g!]), (|41bp . and (|41cp ]. 
The system looks heavily overdetermined. Normally, this 
problem does not arise. If the Zauderer method is applied 
to a typical boundary value problem, for example, solving 
the Laplace equation with Dirichlet boundary conditions 
on part of the boundary and Neumann conditions on the 
remainder [45j , then we will obtain boundary conditions 
for some combinations of variables at the first bound- 
ary and for others at the second, with the total number 
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of conditions just corresponding to the total number of 
equations. However, our problem is not typical, as is 
well-known. The interface position itself is an unknown 
of the problem, requiring the imposition of an additional 
boundary condition at the interface. As a consequence, 
wc obtain a full set of boundary conditions already from 
(|45p . but it is in terms of the unknown interface posi- 
tion /i(^). The remaining boundary conditions p6| then 
are solvability conditions to be imposed on that unknown 
function. It will turn out that three of these conditions 
can be satisfied automatically by requiring h{S.) or the 
curvature to approach zero sufficiently fast at infinity. 
The last one is a non- trivial equation for h{Q which re- 
places the integro-differential equation derivable in prob- 
lems with linear bulk equations. Applying the Kruskal- 
Segur method to this interface equation, we may then 
derive the selection equations. 



VI. SOLUTION OF THE DECOMPOSED 
EQUATIONS 

All equations to be solved are now first order with lin- 
ear derivative terms. This suggests to try their analytic 
solution via the method of characteristics, a step allow- 
ing to make progress despite the nonlinearity of the basic 
equations. 

The principal parts of Eqs. (|4l]) correspond to two 
sets of characteristic coordinates. We start with (|41cp 
and (|44ap . first setting x = x{s,t) with s = s(^,r?) and 
T = t(^, 77). The linear combination of derivatives should 
correspond to a derivative with respect to s only, which 
yields the characteristic equations 



ds ' ds 



^=0. 

ds 



(47) 



Solving this system with the initial condition 

77(5 = 0) = 1, as = o) = T 



we obtain 



-i('7-l) , 
e + i(^/-l) , 

iPf [Th{r)]' 



X = x(t) 



(48) 

(49a) 
(49b) 

l.ik' ' 

i.e., X is simply the analytic continuation into the upper 
77 half plane of the function represented by the boundary 
condition at 77 = 1 (|45dp . This was to be expected, since 
Eq. (|41c|) contains only a principal part. Analyticity re- 
quires X to remain bounded for ^ 00; in fact we have 
the stronger condition (j46dp . To make sure it is satis- 
fied, we may impose the perturbation /i(r) to decay fast 
enough for r ^ ioo so that [Th{T)]' 0. 

The case of the function M is more interesting 
have the same characteristic coordinates s and r 
the equation for M takes the form 



We 
and 



M, 



P. 



(2s + t) M + 



iPf[Th{T)]' 



0, 



2 ^ - ■ 2(1 + 1/ 7') 

Solving this with initial condition (|45ap . we find 



(50) 



A'/(s,t) 



[(l + ir)fe]' 
(l + i/i') 



t{1 + W)\ I 



(51) 



and boundary condition (|46ap is satisfied, if again we 
assume /i(t) -> for r ioo. 

The characteristic coordinates for the other two equa- 
tions are 



5 = 1(77-1) , T = ^-i(7/-l) 
and we obtain 



dN' 
~df 

giving the obvious solution 



0. 



iV" = - 



1 cr[K(f)a(6')]' 

2 l-i;7'(f) 



(52) 



(53) 



(54) 



and boundary condition (j46cp is satisfied, if we require 
the (derivative of the) curvature term na to vanish for 
-ioo [ill. 

Finally, the equation for N becomes 



— fAf(-s,f + 2s) + 



iP/[(f+2s)/7(f + 2s)]'^ii,2 



2(l + i/i'(f + 2s)) 



(55) 



with the boundary condition at s = 0, following from 
(|45bjl and (|45cll with ^ = f : 

^ 0) = - i KizilMl)! . (56) 

^ ' 2 1-i/7'(f) 2 l-ih'{f) ^ ' 

Equation ((55|) can be solved by direct quadrature, with 
the result: 



N{s,f)=-^f J M {-s,f + 2s)ds 



\(f + 2s')h(f + 2s')]' 

e 2 * ds 



LP/ r [{f + 2-s') 
2 J l + ih' 



(f +2s') 



a(«a)'(f) [(l-if)/7(r)]' 



l-i/7'(f) l-ih'{f) 



(57) 



and a sufficient condition for boundary condition (j46bp 
to be satisfied is 



lim 7V(s,f) = 0. 



(58) 



Evaluation of this requirement will produce the central 
equation, to which the Kruskal-Segur method can 
be applied. Note that Eq. ((58|) is an equation for the 
interface position h{f) that has to be satisfied identically 
in the single complex variable f. The next task is to 
cast this equation into a useful form. Since this is purely 
technical, we relegate the calculation to the appendix. 
The resulting interface equation is 
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8 S 



Pf 



(59) 



where we have renamed f into ^ for convenience and re- 
placed the argument 9 of the anisotropy function also by 
^ (the dependence a(^) is given in the appendix). Note 
that we can immediately read off the limit of vanishing 
flow (Pf 0) and verify that it agrees with the local 
equation . [This may of course already be seen from 
Eqs. ((57|) and ([58]) .] The full equation is nonlocal but it 
is tractable via asymptotic methods. 

In principle, Zauderer's scheme may be used to solve 
the system of equations ([57]) perturbatively. To carry out 
this (complicated) calculation, one would have to keep 
the dropped terms of order and to add terms contain- 
ing the missing eigenvectors (and a factor e) to Eqs. (j39bp 
and p9cp . Inspection immediately reveals that this ex- 
pansion in powers of e does not correspond to a pertur- 
bation series about the flow-Ivantsov solution: setting 
e = does not give us the full flow-Ivantsov solution but 
only the solutions of the first-order equations obtained 
from the projection onto eigenvectors of A; e.g., in the 
case of the Laplace equation for x would just obtain 
a solution to X(, + iXr? = 0. However, these pieces of the 
full solution are the ones that diverge in the limit cr ^ 
near the singularity of interest, whereas the other terms 
remain finite. Hence, the lowest-order Zauderer solution 
corresponds to the exact solution of the problem near the 
singularity. So the perturbative scheme arising from Za- 
uderer decomposition corresponds to an expansion about 
the analytic continuation of the fiow-Ivantsov solution 
in the vicinity of the singularity. This may be seen as 
the deeper reason why a condition for the transcendental 
term which is beyond all orders in regular perturbation 
theory appears already at first order in our approach ~ 
near the singularity this term is not small, so it has to 
be present in a Zauderer type perturbation theory. 

We have carried this perturbative approach beyond 
first order for the simpler problem without flow. If we ex- 
pand h and k in powers of e as well, i.e., h = hQ + ehi + . . . 
and K = kq + • • - i a- solvability condition similar to (j34p 
appears to turn up at the next order involving hi and ki 
- but it is automatically satisfied. Hence, it seems that 
the lowest-order solvability condition does indeed capture 
the mismatch function needed to obtain the selection cri- 



terion. While we knew this to be true from comparison 
with known results in the case of ([34|) . these arguments 
suggest it to hold in general, i.e., also for Eq. ([59]) . 

In order to obtain the mismatch function (or the con- 
tribution to it by the singularity considered) at the inter- 
face, we have to solve Eq. ([5^ . Far from the singularity, 
this can be done by linearization in terms of the interface 
position h and its derivatives, which may all be consid- 
ered small. The appropriate tool is WKB analysis. Due 
to the linearity of the problem, this will provide the so- 
lution up to a constant factor only. Using asymptotic 
analysis, we can then solve Eq. (|59p near the singularity, 
taking all important nonlinearities into account. Asymp- 
totic matching of the two solutions provides the prefactor 
of the near-interface solution. The mismatch function 
calculated from it must vanish at the tip of the needle 
crystal - this is the solvability condition. 



VII. WKB ANALYSIS FAR FROM THE 
SINGULARITY 



Linearizing (j59[) . we obtain the inhomogeneous linear 
equation 

/ 1 h"{0 ih'iO ■ 



Pf L 



= (l-iOMe) 

C'Ci-iOMOde'. (60) 



The solution of this consists of a particular solution to 
the inhomogeneous equation (which will be captured by 
regular perturbation theory) plus the general solution of 
the homogeneous equation (with integration constants 
to be determined from boundary conditions on h). The 
latter consists of an exponentially small and an expo- 
nentially large term. The large term is suppressed al- 
ready within regular perturbation theory, but the small 
one will not appear therein at any finite order. It be- 
comes important, when, due to symmetries of the prob- 
lem, all terms of regular perturbation theory vanish. In 
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the needle-crystal problem, this is the case at the tip of 
the crystal. So the transcendentally small term that must 
be suppressed can be identified with the decaying solu- 
tion of the homogeneous linear equation corresponding to 
(|5(H) . Alternatively, we could argue that the general so- 
lution to the inhomogeneous equation may be obtained, 
within WKB theory, via the method of variation of con- 
stants [4l|. Again, the exponentially large term must 
be eliminated by an appropriate choice of an integration 
parameter. The exponentially small one has the same 
form as the decaying solution of the homogeneous equa- 
tion, except that there is now a slowly varying prefactor. 
Since the mismatch function is to be evaluated at the tip 
position in the end, it has the same form as this solution. 

The only tricky part of the calculation of the WKB 
solution is the evaluation of the integral in ([SH]) . which 
can be done via integration by parts. We obtain 



/i(C) = Bie-(l + ie)"Ml-ie)" 
with an unknown constant Bi and 



SqU) 

e 



16 ? 



(61) 



^o(e)=i /(l + iOMl-iO^dC. (62) 



VIII. SOLUTION NEAR THE SINGULARITY 



The most appropriate form of Eq. ([59| for a local anal- 
ysis near the singularity seems to be Eq. (|B9p . With M 
given explicitly, it reads 



^'e*(«'^-«"^)([(i+ir)Mr)]' 



|;[W)]'-(l-c^«"(^^"-«'))) 



de"d^' 



(63) 



Introducing the stretching transformation[47j 

e = -i(l-^"0 (64) 

with a = y (obtained from a dominant balance consider- 
ation), we set 



from which we get 



h' = 

h = -a 

.2q , 



and find, to leading order in a 



(65) 

(66a) 
(66b) 
(66c) 



1 



1 + , 



(2t + 2<?!)) ^ 



(t)t 



(67) 



Rewriting the right-hand side of Eq. ([S^]) in terms of the 
new variable is a bit more involved and leads to 
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+ 2{t' - 1") , 



i+0(t") 

At" At' . 



e 1 



■cr°(t'-t") 



(68) 



An important result is that the leading order of a is cr^" 
which in the limit <t ^ vanishes faster than the lead- 
ing order of Eq. (jBT)) . This means that for finite Pf the 
selected stability parameter will be the same as in the 
flowless case with otherwise equal parameters. Hence, for 
the same Pc, the same tip radius (and velocity) will be 
selected with and without flow. For given undercooling, 
the selected velocity and tip radius will be different from 
the corresponding quantity of the flowless case only due 
to the different relationship (|27p between undercooling 
and growth Peclet number, which contains a dependency 
on Pf (i.e., for the same A, Pc is different in the two 
cases). In fact, this result has been used as an assump- 
tion in the past to compute selected growth velocities in 
convective situations [48|. Here, it has been proved for 
the case of potential flow, but our experience with other 
flow patterns suggests this to be a general feature of con- 
vection. To our knowledge, no general proof has been 
given so far. 

To obtain a nontrivial dependency of the stability pa- 
rameter on the flow Peclet number, we have to assume 
large flow velocities, e.g. Pf = 0((t"") = 0(cr"^/''). Hence 
we set 



= fa". 



(69) 



Since we expect that there is no solution in the 
isotropic case, we take surface tension anisotropy into 
account right away. Performing the stretching transfor- 
mation for a{9), we find 



a{9) = 1 



2/3cr-2" (i 



{t + <pf{l + ^y 



(70) 



The anisotropy parameter /3 usually is numerically small. 
Setting 

f] = a^^b = a^'^b, (71) 
we may cast our interface equation into the form 



<t)t + Pi 



1 + ^{t") 



,Piit'-t") 
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(72) 



Given the boundary condition that the imaginary part 
of (j) vanishes for t ->• oo (which is the condition that 
the tip slope of the needle crystal is equal to zero) and a 
prescribed value of Pi , this constitutes a nonlinear eigen- 
value problem for b. We have solved this numerically in 
the complex plane, using a scheme similar to the one 
given by Tanvccr [23|; we employ a relaxation method 
along two straight intersecting lines in the complex plane, 
one of them parallel to the imaginary, the other lying on 
the real t axis. Details of the numerical approach, which 
is a root finding problem involving the integration of sev- 
eral ordinary differential equations and exhibits a certain 
level of complexity, will be given elsewhere. 

We do not find any solutions with 6 = 0, suggesting 
that there does not exist, as anticipated, any steady-state 
needle crystal close to a flow-Ivantsov parabola in the 
case of isotropic surface tension. 

For anisotropic surface tension, we have the usual re- 
lationship between the selected stability parameter and 
the anisotropy parameter 



(73) 



If the general solution behavior is similar to that of the 
flowless case, the solution corresponding to the lowest 
eigenvalue b should be the only one that is linearly stable. 
We assume this to be true, but have not yet been able to 
check it. 

The relationship between the physical flow Peclet num- 
ber and our Pi is 



Pf = 4Pn 



(74) 



Once we have Pf and cr, we may determine Pc (numer- 
ically Hi]) from Eq. 1^ and using the definitions ^ 
of Pc and (pn)) of cr we can evaluate both the selected tip 
radius p and tip velocity V. 

Note that while our approximations hold in the limit 
Pc « Pf, which implies in particular an approximation 
for T' in Eq. ([26|) that does not approach the limit Pf 
uniformly in 77, the eigenvalue b obtained numerically 
will still be correct in that limit, due to the structure of 
Eq. (|59p which reduces to the selection criterion without 
flow. Indeed, we have verified that we obtain the same 
value of b as Tanveer [2^ in the case without flow. 

Although our model is definitely a toy model - ex- 
perimental flow patterns and velocities will not be well 
described by a potential flow [H^ - we carry the cal- 
culation to its end using parameters determined for an 
experimental substance, pivalic acid. Since it is not to be 
expected that this will give more than qualitative trends. 



the purpose of this exercise is mostly to demonstrate that 
the (relatively elaborate) formalism produces numbers fi- 
nally and that these numbers do not have unreasonable 
orders of magnitude. 

Caveats to be kept in mind are: 
- We use the symmetric model, whereas the one-sided 
model would be more appropriate for experiments with 
solute diffusion. However, this is known to just make a 
difference of a factor of two in the selected velocity [5l[ in 



the diffusion-limited case. We expect a similar closeness 
of results of the two models in the presence of convection. 

- Our model is only two-dimensional, which certainly im- 
pedes its quantitative applicability to experiments. On 
the other hand, typically the predictions of microscopic 
solvability theory do not differ much for two-dimensional 
and (axisymmctric) three-dimensional systems [s^ . 

- More importantly, pivalic acid has kinetic anisotropy, so 
it is not to be expected anyway that a model imposing lo- 
cal equilibrium at the interface will yield a good descrip- 
tion. We chose the experiments from Ref. [s^ for com- 
parison, because they have flow velocities that are in the 
range of numerical accessibility for our code, whereas in 
experiments with succinonitrile j54| (a system expected 
to be better suited for comparison on physical grounds), 
the imposed flow velocities were very large, leading to 
convergence problems in our eigenvalue computation. 

- Potential flow and hence our relationship between Pc 
and Pf is not realized in the experiments. 

Material parameters were taken from Refs.[5^[5^ and 
[56I and an undercooling of about 0.2 K (equivalent to 
A = 0.0169) was assumed, corresponding to a situation 
considered in the experiments. Results are shown in 
Figs. [1] to HI 
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FIG. 1. The stability parameter a as a function of the flow 
Peclet number Pf. Material parameters used correspond to 
pivalic acid with A = 0.0169 and 13 = 0.08 HI]. 

Figure [T] gives the selected value of cr as a function 
of the flow Peclet number for fixed undercooling A and 
an anisotropy parameter that corresponds to a measured 
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value [5£ 
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FIG. 2. The crystal growth velocity as a function of the 
flow velocity U. 
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FIG. 3. The tip radius of the needle crystal p as a function 
of the flow velocity U. 



Finally, we do compare the values of p^V obtained from 
this calculation with experimental values in a flow situ- 
ation [5^, because p^V is expected to be a slowly vary- 
ing quantity and therefore what matters mostly is the 
overall order of magnitude. As Fig. 2] demonstrates, this 
quantity compares reasonably with experiment. In fact, 
considering that the experimentalists describe their flow 
pattern as approximate Oseen flow, the agreement is not 
too bad. This should of course not be taken too seriously 
either. A real comparison will have to await a calcula- 
tion with a more realistic flow (and, for pivalic acid, a 
different interface boundary condition). 

To conclude this section, it may be noted that a lo- 
cal asymptotic analysis of Eq. ([72|) for t ^ 00 yields the 
same transcendental behavior as Eq. (|61|) and provides 
the prefactor Bi in terms of the solution of the nonlin- 
ear equation (|72p as a function of Pi and b. Since the 
boundary condition on the tip slope was however already 
incorporated into the numerical scheme for the solution 
of Eq. ([7^ . this calculation does not provide anything 
new. 
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FIG. 4. The product as a function of the flow velocity 
U, compared with an experiment [ssj . 



In Figs. [2] and [3] we give the selected growth velocity 
V and tip radius p [s^ in dimensional form. We refrain 
from comparing these numerical results with a concrete 
experiment, because there are too many uncertainties re- 
garding the applicability of the toy model to real life. 

All that we wish to point out here is that there are 
power law scaling relations between the growth velocity 
and the velocity of the imposed flow as well as between 
the tip radius and the flow velocity, valid in a range of un- 
dercoolings. This feature will probably not disappear in 
a more quantitative calculation. In fact, we have checked 
for an extended range of anisotropies, thus varying a be- 
tween small and very large values, that the scaling expo- 
nents change only slightly. 



IX. CONCLUSIONS 

After introducing the combination of Zauderer de- 
composition with the Kruskal-Scgur approach recently 
[3^ Issj , we have now presented the method in more de- 
tail. The analytic part of the calculation has been exem- 
plified with a fully nonlinear problem. Approximations 
that were introduced in [s^] for didactic reasons have 
been removed, rendering the full power of the method 
visible. 

We believe our approach to be the only one presented 
so far that has the potential of a rigorous solution of pat- 
tern selection problems with nonlinear bulk equations. 
Essentially, our belief that the method is rigorous rests on 
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two facts: first, the Zauderer decomposition scheme pro- 
duces a solution that becomes exact near the appropriate 
complex-plane singularity; second, the WKB solution de- 
rived from the interface equation within the scheme gen- 
erates the same transcendental terms that a WKB solu- 
tion derived from the full system of partial differential 
equations would. The second statement has been shown 
to be true for the flowless case [11] and we have given ar- 
guments here, why it should carry over to the nonlinear 
case as well. 

The elegance and power of the method show up in its 
rendering the purely diffusion-limited case almost triv- 
ial [3^ Issj . When applied to a problem with nonlinear 
bulk equations, calculations certainly become involved. 
But the problem remains solvable in a controlled man- 
ner, not provided by other methods. That in the final 
step the numerical determination of an eigenvalue be- 
comes necessary should not prevent us from considering 
the approach basically analytical. A similar final step 
is necessary in almost all related problems with simpler 
bulk equations, even though the nonlinear equation to be 
solved numerically is less difhcult in these cases. 

We are convinced that our method will render a num- 
ber of selection problems accessible to solvability theory 
for which controlled approximations could not be devel- 
oped in the past, thus opening a new line of research. 
These problems would include nonlinear diffusion jssj . 
density-driven convection [soj (for which we ha ve g iven 
a preliminary treatment before [6d|). Oseen flow |53l. [6l|. 
the effect of the Kapitza resistance on dendritic growth 
of helium (62l . |63| . but also pattern selection problems 
outside of crystal growth such as, for example, the mo- 
tion of the two-phase front between superconducting and 
normal conducting parts of a material 64 1 . 
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Appendix A: Conformal parabolic coordinates 



The unnormalized coordinate basis is given by 



ax _ (^) 



(Al) 



which yields 5| e \£^f = e + = l^vf = gfy 



For the nabla operator, we get 



V 



1 



d d 



whereas the Laplacian reads 
1 



(A2) 



^2 _,_ ^2 



(A3) 




= 3 



FIG. 5. Coordinate lines for conformal parabolic coordinates. 
Note that these coordinates will cover the xy plane twice, if 
negative values for rj are admitted. 



After normalization, the basis vectors are 
1 



(77e^, -^e.y) , 



(A4) 



which can be inverted to express the Cartesian basis by 
the orthonormal parabolic one 



1 



(jye^ +^e^) 



(A5) 



Let 'ris{£,) - rj = Q describe the interface, then the normal 
vector n can be derived from Frenet's formulas. The 
position vector at the interface may be written 



1 



(A6) 



and the differential line element along this curve is 

ds = Vdx^ + d?/2 = ^/{r]1+£,^^){l + r]'^^)A£,. (A7) 

The tangential vector at the interface is given by dx/ds, 
the normal vector must be orthogonal to it. By this con- 
dition, it is determined up to a sign that we choose so 
as to make the normal vector point into the liquid. This 
procedure yields 



n : 



(A8) 
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The curvature is given by 
d^r 



ds2 
1 



Vs 



(A9) 



and it is positive for a convex solid. 

We assume the usual model of four-fold crystalline 
anisotropy: 

a(6i) = 1 - ^ cos 46* = 1 - /? ( 1 - 8 cos^ 6 sin^ 9) (AlO) 

The small parameter /3 is the strength of the anisotropy. 
is the angle of the interface normal with the y axis, 
so we have cos 9 = n - ey and sin 9 = n • e^; , which allows 
us to find the anisotropy function expressed in parabolic 
coordinates. 



i{9) = 1-/3 



1- 



(All) 



Finally, the flow velocity is given by 
1 



Appendix B: Derivation of the interface equation 



We first introduce some simplifications of notation. 



Substituting t = ^ and s' = ^(u - ^) and defining 



z{0 = i-ih'{0, 



(Bl) 
(B2) 
(B3) 



we have from Eq. f|57|) with f|58|l 

ioo 

F'io = -^mm /A/(i(^-w),w)d« 

ioo / 

2 J z{u) 



(B4) 



with the prime denoting a derivative with respect to ^ or 
M, depending on whether the term concerned is outside 
or inside an integral. Writing out M, we have 



M 



2z{u) 
V 2 [uh]' 



[(1 + iu) h]' - =-^^^ (l - e-^(«-")' 
u \ 



(B5) 



and 



_!^Me*«-)\ (B6) 
4 z{uj 

Inserting this into (jB4p . we obtain a useful expression for 
the derivative of F: 

F'(0=-2iz-(0 /^M(i(^-«),«)du 



-2iz{0 



loo 



IOO 



This can be integrated by parts. Using 



-2i— 



(B7) 



(B8) 



we arrive at 



IOO 

(A12) F{0=-2iz{0 f Mi^^{^-u),uyu 



^(0 



^(0 



(B9) 



which is not quite the form we want. On the one hand, 
equation (|B9[) manifests a certain generality, since it is 
valid for arbitrary flows. But on the other hand, we 
would appreciate to have a right hand side that obvi- 
ously vanishes for Pj ^ 0. To achieve this, we use (jB9P 
to eliminate the first term on the right hand side of (|B4p : 



4 4 z{0 



2 

Pf 



t, loo 

e / ^"iO f ^/(i(e'-^i),^i)d^ide' 

^(^) r i!l!^e*(«-")'dii (BIO) 
J z(u) 



Employing the identities 



loo 

/ 



z{u) 



du ■ 



^(0 
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-/ 



du . 



F{Oe- 



(Bll) 
, (B12) 



the first of whicli is obtained via integration by parts 
again, we may rewrite (|B10p as follows 



iz{u) 



) 



du 



El 

2 



e / h"{a f M(i(e'-u),u)dud^'. 



(B13) 



d 

e 8 « — 



4 z{0 



With one further integration, we arrive at Eq. ((59 
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